function [B] = generatespy_index(Nx,Ny,  m_ccc,m_css,m_cww,m_cee,m_cnn)

%Generated the Pentadiagonal Matrix of size Nx*Ny x Nx*Ny
%   Detailed explanation goes here

B = zeros(Nx*Ny, Nx*Ny);

for ii = 2:Nx+1
    for jj = 2:Ny+1
        index = (jj-2)*Nx+(ii-1);
        B(index,index) = m_ccc(jj,ii);
    end
end

%% -- Inner Points ------
for ii = 3:Nx
    for jj = 3:Ny
        index = (jj-2)*Nx+(ii-1);
        B(index,index-1)= m_cww(jj,ii);
        B(index,index+1) = m_cee(jj,ii);
        B(index,index+Nx) = m_cnn(jj,ii);
        B(index,index-Nx) = m_css(jj,ii);
    end
end
%% North and South Boundaries
for ii = 3:Nx
    % -- South Boundary
    jj=2;
    index = (jj-2)*Nx+(ii-1);
    B(index,index-1)= m_cww(jj,ii);
    B(index,index+1) = m_cee(jj,ii);
    B(index,index+Nx) = m_cnn(jj,ii);
    % -- North Boundary
    jj=Ny+1;
    index = (jj-2)*Nx+(ii-1);
    B(index,index-1)= m_cww(jj,ii);
    B(index,index+1) = m_cee(jj,ii);
    B(index,index-Nx) = m_css(jj,ii);
end
%% -- East and West Boundary
for jj = 3:Ny
    ii=2;
    % -- West Boundary
    index = (jj-2)*Nx+(ii-1);
    B(index,index+1) = m_cee(jj,ii);
    B(index,index+Nx) = m_cnn(jj,ii);
    B(index,index-Nx) = m_css(jj,ii);
    % - East Boundary
    ii=Nx+1;
    index = (jj-2)*Nx+(ii-1);
    B(index,index-1)= m_cww(jj,ii);
    B(index,index+Nx) = m_cnn(jj,ii);
    B(index,index-Nx) = m_css(jj,ii);
end


%% -------- Corner Points ---------%%
% South West Corner Point
ii=2;
jj=2;
index = (jj-2)*Nx+(ii-1);
B(index,index+1) = m_cee(jj,ii);
B(index,index+Nx) = m_cnn(jj,ii);

% North West Corner Point
ii=2;
jj=Ny+1;
index = (jj-2)*Nx+(ii-1);
B(index,index+1) = m_cee(jj,ii);
B(index,index-Nx) = m_css(jj,ii);

% North East Boundary
ii = Nx+1;
jj = Ny+1;
index = (jj-2)*Nx+(ii-1);
B(index,index-Nx) = m_css(jj,ii);
B(index,index-1)= m_cww(jj,ii);

% South East Boundary
ii=Nx+1;
jj=2;
index = (jj-2)*Nx+(ii-1);
B(index,index+Nx) = m_cnn(jj,ii);
B(index,index-1)= m_cww(jj,ii);



%% -----------------



%
% A = zeros(Nx3*Ny3,Nx3*Ny3);
% for ii=2:Nx+1
%     for jj= 2:Ny+1
%         index = (jj-1)*(Nx3) + ii;
%         A(index,index) = m_ccc(jj,ii);
%         A(index,index-1)= m_cww(jj,ii);
%         A(index,index+1) = m_cee(jj,ii);
%         A(index,index+Nx3) = m_cnn(jj,ii);
%         A(index,index-Nx3) = m_css(jj,ii);
%     end
% end
%
% for ii=2:Nx+1
%     for jj= 2:Ny+1
%         index = (jj-2)*Nx+(ii-1);
%         global_index = (jj-1)*(Nx3) + ii;
%         B(index,index) = A(global_index,global_index);
%         if (index >1 )
%
%         B(index,index-1) = A(global_index,global_index-1);
%         end
%
%         if (index <(Nx*Ny-Nx))
%         B(index,index+1) = A(global_index,global_index+1);
%         end
%
%         if (index <(Nx*Ny-Nx))
%
%         B(index,index+Nx) = A(global_index,global_index+Nx3);
%         end
%
%         if (index > Nx)
%         B(index,index-Nx) = A(global_index,global_index-Nx3);
%         end
%
%
%     end
% end
%
%
%
%
% % %spy(A);
% % k = Nx3+2;
% % l = Ny*Nx3+(Nx+1);
% % B = A(k:l , k:l);
% % %B = A(4:(Nx+1)*(Ny+1),4:(Nx+1)*(Ny+1));


end

